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ABSTRACT 


By  using  a  linear  systems  theory  approach,  an  ocean 
medium  transfer  function  based  on  the  WKB  approximation  can  be 
developed.  The  phase  computations  for  the  transfer  function 
are  made  by  evaluating  the  WKB  phase  integral . 

Two  applications  of  ray  acoustics  theory  are  investigated 
as  accurate,  efficient  alternatives  to  direct  numerical 
integration  of  the  WKB  phase  integral.  Both  applications  base 
phase  computations  on  signal  travel  time.  The  difference  is 
their  treatment  of  the  sound-speed  versus  depth  data  pairs. 
One  forms  a  sound-speed  profile  by  using  the  piecewise  linear 
approximation  method  while  the  other  uses  an  Akima  cubic 
spline  fit  to  the  data. 

Each  method  can  identify  source-to-receiver  eigenrays  and 
provide  ray  trace  plots. 
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I .  INTRODUCTION 


A.  BACKGROUND 

In  analyzing  ocean  acoustic  pulse-propagation  problems, 
accurate  phase  calculations  roust  be  performed  in  order  to 
predict  the  acoustic  signal  at  the  receiver.  Linear  systems 
theory  provides  analytical  expressions  for  analyzing  the 
propagating  acoustic  field  [Refs.  1-4].  These  well-known 
equations  form  the  basis  for  the  FORTRAN  programs  developed  in 
this  thesis. 

Program  input  is  depth  versus  sound-speed  data  pairs  where 
the  speed  of  sound  is  an  arbitrary  function  of  depth.  The 
ocean  is  viewed  as  a  linear,  time-variant,  space-variant 
filter.  The  WKB  approximation  can  specify  this  filter's 
transfer  function  [Refs.  1-4].  For  the  arbitrary  sound- 
speed  profile,  no  exact  transfer  function  exists.  The 
transfer  function  requires  a  method  for  evaluating  signal 
phase  at  the  receiver. 

Two  phase  evaluation  methods  are  presented  and  contrasted. 
Each  method  calculates  the  phase  of  the  acoustic  signal  for  a 
specified  horizontal  range.  The  first  method  overlays  a 
piecewise  linear  profile  on  the  input  sample  values  for  the 
speed  of  sound  at  various  depths.  The  theory  of  ray  acoustics 
is  used  to  calculate  travel  time  and  phase.  The  second  method 
fits  a  smooth  cubic  spline  curve  to  the  input  samples.  This 
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method  solves  the  propagation  problem  using  a  system  of  three, 
first-order  differential  equations  [Ref.  5],  Both  methods 
make  phase  calculations  using  signal  travel  times. 
Additionally,  each  method  is  capable  of  identifying  eigenrays 
or  rays  that  directly  connect  the  signal  source  to  the 
receiver. 

Another  method  is  introduced  to  validate  results.  Direct 
numerical  integration  is  performed  to  calculate  phase  for  a 
single  gradient,  linear,  sound-speed  profile.  A  separate 
FORTRAN  program  implements  this  numerical  integration  routine 
providing  a  totally  independent  verification. 

An  overview  of  the  theory  behind  each  method  is  presented 
in  Chapter  II.  Computer  simulation  results  are  presented  in 
Chapter  III  for  various  input  sound-speed  profiles.  Tabular 
results  of  phase  calculations  are  presented  for  each  of  the 
three  analysis  methods.  Using  input  from  the  piecewise  linear 
and  the  cubic  spline/differential  equation  solution  methods, 
ray  traces  of  the  propagating  field  are  presented  as  a  visual 
aid  in  interpreting  the  results. 
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II.  THEORETICAL  BACKGROUND  FOR  THE 
EVALUATION  OF  PHASE  INTEGRALS 

A.  UNDERWATER  ACOUSTIC  PULSE  PROPAGATION 

Ocean  acoustic  pulse-propagation  models  can  be  derived  by 
using  the  principles  of  linear,  time-variant,  space-variant, 
systems  theory  and  the  physics  of  wave  propagation  in 
inhomogeneous  media.  Linear  systems  theory  allows  for  the 
development  of  an  ocean  medium  transfer  function.  An  ocean 
medium  transfer  function  that  is  based  on  the  WKB 
approximation  has  been  derived  and  is  given  by  References  1 
and  2  as  follows: 


H  (f ,  fr/Yo'Y) 


A 

\/\kyiy)  ■ 


e 


-)  f  t<,(C)ilC 
J»o 


gJ2'Tf,(y-y„) 


(2.1) 


where 


f  is  the  frequency  in  Hz 

f^  is  the  radial,  spatial  frequency  in 
cycles  per  meter 

yg  is  the  source  depth  in  meters 


A 


kyiYo) 


(2.2) 
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ky(y)  is  the  propagation  vector  component  in 
the  Y  direction  with  units  of  radians 
per  meter  and  is  given  by 


■±  2n{[f/c{y)f  -  <  f/c(y) 

j27r{f,^  -  [f/c(y)]^}'^  ,f,  >  f/c(y) 


(2.3) 

(2.4) 


c(y)  is  the  sound-speed  expressed  as  a  function 
of  depth  with  units  of  meters  per  second 


fy  is  the  transmitted  (input)  spatial  frequency 
in  the  Y  direction  at  the  source  in  cycles  per 
meter  and  is  given  by 

r±  [(f/Co)^  -  ,  f^  <  f/Co  (2.5) 

[t  [f^2  -  (f/Co)^]"'  ,  f^  >  f/Co  (2.6) 

and 


Cq  is  the  speed  of  sound  in  meters  per  second  at 
the  source  depth  yQ,  that  is,  Cq  =  c(yQ); 

The  plus  (minus)  sign  in  Eqs.  (2.3)  and  (2.5)  is  chosen 
whenever  y  -  yo  >  0  (y  -  Yq  <  0)  •  The  minus  (plus)  sign  in 
Eqs.  (2.4)  and  (2.6)  corresponds  to  the  plus  (minus)  sign  in 
Eqs.  (2.3)  and  (2.5) .  This  thesis  deals  only  with  propagating 
waves.  Therefore,  Eqs.  (2.4)  and  (2.6)  representing  the 
generation  of  evanescent  waves  (i.e.,  decaying  exponentials), 
will  not  be  used. 

The  evaluation  of  the  phase  integral  in  Eq.  (2.1),  namely. 


ey(y) 


(2.7) 
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has  been  attempted  in  past  studies  in  which  solution 
techniques  included  direct  integration  and  binomial  expansions 
[Ref.  3]. 

Computer  programs  using  these  techniques  proved  to  be  very 
expensive  in  terms  of  CPU  time.  This  thesis  compares  four 
different  methods  for  evaluating  the  phase  integral.  The  four 
methods  are  presented  in  the  following  sections. 

B.  DIRECT  INTEGRATION 

Computer  code  was  developed  to  directly  evaluate  the  phase 
integral  given  by  (Eg.  2.7)  using  numerical  integration 
routines  from  the  International  Math  and  Statistics  Library 
(IMSL) .  Direct  integration  is  a  time-consuming  technique  that 
was  applied  only  under  the  following  constraints: 

-  free-space  acoustic  propagation 

-  linear  sound-speed  profile  with  a  single  gradient  g 
in  inverse  seconds  and 

-  propagating  signals  that  have  not  passed  through  a 
turning  point  in  the  medium. 

Expressions  for  the  propagation  vector  component  ky(y)  of 
the  phase  integral  are  given  by  Eqs.  (2.3)  and  (2.4).  For  a 
linear  sound-speed  profile  with  a  single  gradient,  the 
dependence  of  c(y)  on  depth  y  is  given  by 

c(y)  =  c(yo)  +  g(y  -  y^)  (2.8) 

where 

c(yQ)  is  the  sound-speed  at  the  source  depth 
yQ  meters 
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g  is  the  sound-speed  gradient,  and 
y  is  the  desired  depth  for  sound-speed. 

Since  the  direct  integration  approach  is  limited  to  free- 
space,  single  gradient,  propagation  problems,  it  is  used  only 
to  validate  the  results  of  other  solution  techniques  for 
several  simple  test  cases.  The  results  will  be  compared 
against  the  next  two  methods  to  be  presented.  These  methods 
will  apply  the  theory  of  ray  acoustics  [Ref.  4  &  5]  to  find 
travel  time.  Phase  calculations  will  be  based  on  travel  time 
calculations . 

C.  PHASE  CALCULATIONS  BASED  ON  PIECEWISE  LINEAR  SOUND-SPEED 

PROFILES 

A  FORTRAN  program  developed  by  Lim  [Ref.  5]  was  the  ray 
acoustics  propagation  code  that  was  used  to  draw  ray  trace 
plots  and  to  calculate  travel  time,  in  addition  to  the  depth, 
angle  of  propagation  and  path  length  along  a  ray  path  as  a 
function  of  horizontal  range.  This  method  applies  ray 
acoustics  to  a  piecewise  linear  sound-speed  model  of  the  ocean 
medium.  Based  on  the  input  sound-speed  versus  depth  data 
pairs,  layers  are  defined  in  the  ocean  medium.  Each  layer  has 
an  upper  and  lower  boundary  at  specific  depths.  A  constant 
sound-speed  gradient  g  is  calculated  for  each  layer.  The 
sound-speed  for  any  desired  depth  is  computed  using  Eg.  (2.8) 
with  the  appropriate  gradient. 

In  a  constant  gradient  medium,  ray  acoustics  theory  allows 
calculation  of  travel  time  with  closed-form  equations.  For  an 
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incremental  increase  in  horizontal  range  from  the  source,  the 
angle  of  arrival  for  the  propagating  ray  can  be  shown  to  be 
[Ref.  5]. 

;0(y)  =  cos  ^  [cos;8q  -  b*g(y)  *  rngstp]  (2.9) 

where 

/3q  is  the  ray  launch  angle, 

g(y)  is  the  gradient  in  the  layer  at  depth  y,  and 

rngstp  is  the  incremental  increase  in  horizontal 
range . 

The  ray  parameter  b  in  Eq.  (2.9)  is  given  by  [Refs.  4  &  5]. 

b  =  sin/3o  /  c(yo)  (2.10) 

Knowing  the  arrival  angle  /3(y),  the  ray  depth  y  is  given  by 
[Refs.  4  &  5] 

^0 

y  =  Yo^— 
g 


sin/3  (y) 
sin^(yj 


(2.11) 


The  arrival  angle  also  allows  calculation  of  travel  time  from 
[Refs.  4  &  5] 


tan[/3(y)/2] 

tan/3(yQ)/2 


(2.12) 


Once  the  travel  time  (in  seconds)  is  found  for  the  desired 
horizontal  range,  total  phase  is  given  by 

0^  =  27rfr  (2.13) 

The  phase  integral  given  by  Eq.  (2.7)  represents  a  phase 

change  in  the  depth  (y)  direction.  The  constant  value  of  the 
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propagation  vector  component  in  the  radial  direction  allows 
for  easy  calculation  of  phase  change  in  the  radial  direcrion. 
This  value  is  given  by 

©„  =  27rf^  *  hrzrng  (2.14) 

where  hrzrng  is  the  total  horizontal  range  traveled. 

The  phase  in  the  depth  direction  is 

Gy  -  ©T  -  ©R.  (2.15) 

Phase  expressed  in  radians  is  a  modulo  In  function.  The 
proper  solution  for  the  phase  integral  Eg.  (2.15)  must  be 
expressed  as  a  modulo  27r  function: 

©y  =  modulo  (©y  ,  2n  )  .  (2.16) 

This  method  (referred  to  as  method  1  in  Ref.  5)  is  capable 
of  analyzing  propagating  rays  at  any  horizontal  range  in 
either  free  space  or  a  bounded  medium. 

D.  PHASE  CALCULATIONS  BASED  ON  AKIMA  CUBIC  SPLINES  AND 
ORDINARY  DIFFERENTIAL  EQUATIONS 

This  application  of  ray  acoustics  applies  a  more 
sophisticated  treatment  to  the  depth  versus  sound-speed  data 
pairs  sampled  from  the  ocean  medium.  The  data  pairs  are  used 
to  form  a  smooth  Akima  cubic  spline  sound-speed  profile. 
Splines  offer  the  advantage  of  using  all  data  points  in 
generating  a  profile,  and  they  place  no  restrictions  on 
spacing  between  data  points.  The  Akima  version  of  the  cubic 
spline  was  chosen  for  its  excellent  ability  to  combat  wiggles 
in  the  profile,  that  is,  it  suppresses  oscillations  that  would 
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cause  overshoots  and  undershoots  in  the  sound-speed  versus 
depth  profile. 

This  method  uses  ray  acoustics  theory  to  generate  a  system 
of  three,  first-order  differential  equations  [Ref.  5]: 

Yi  =  72  (2.17) 


(2.18) 


^3 - i -  (2.19) 

(y^) 

where 

c(y^)  is  the  derivative  of  the  sound-speed  with 
respect  to  depth  at  depth  y, , 

y,  is  the  ray  depth, 

yj  is  the  cotangent  of  the  ray's  arrival  angle 
/3(y)  ,  and 

yj  is  the  travel  time  of  the  ray. 

Once  solved,  this  system  of  differential  equations  allows 
phase  to  be  calculated  from  travel  time.  The  phase 
calculation  is  performed  exactly  as  shown  for  the  piecewise 
linear  sound-speed  profile. 

Like  the  previous  method,  the  Akima  cubic 
spline/differential  equation  solution  can  be  used  on  free 
space  and  bounded  media  problems. 
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E.  EIGENRAYS 


Eigenrays  are  propagating  rays  that  exactly  connect  the 
sound  source  to  the  receiver.  The  FORTRAN  propagation  code 
developed  can  search  for  and  identify  eigenrays.  The  input 
required  is 

-  depth  versus  sound-speed  data  pairs, 

-  source  depth  in  meters, 

-  receiver  depth  and  range  in  meters, 

-  angle  step  in  degrees  between  possible  eigenrays  to 
be  evaluated,  and 

-  allowed  depth  error  denoted  ye^ror  in  meters. 

Rays  passing  within  the  allowed  error  or  tolerance  ygrror 
receiver  are  identified  as  eigenrays.  The  eigenray  mode  can 
employ  either  the  piecewise  linear  sound-speed  profile  or  the 
Akima  cubic  spline/differential  equation  method  for  ray 
propagation. 

The  eigenrays  are  found  by  trial-and-error .  This  search 
method  was  chosen  after  the  IMSL  DBVPMS  program  failed  to 
solve  the  problem.  The  DBVPMS  program  is  a  differential 
equation  solver  that  was  applied  to  the  system  of  differential 
equations  given  by  Eqs.  (2.17)  through  (2.19).  The  routine 
uses  the  shooting  method  to  find  eigenray  solutions  to 
boundary  value  problems.  It  was  unable  to  converge  to  a 
solution.  Solutions  for  the  acoustic  problems  investigated 
are  difficult  for  this  algorithm  because  of 


10 


-  long  propagation  distances 


-  the  inhomogeneous  ocean  presents  a  continuously 
varying  medium  with  discrete  boundaries,  and 

-  the  system  of  differential  equations  has  non-constant 
coefficients  whose  values  are  a  function  of  depth. 

When  the  eigenray  mode  is  selected,  rays  that  are 

refracted  or  reflected  to  within  y^pror  of  the  receiver 

discarded. 


not 

are 
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III.  COMPUTER  SIMULATION  RESULTS 


A.  OVERVIEW 

The  computer  simulation  results  presented  in  this 
chapter  perform  phase  calculations  intended  to: 

-  show  that  ray  acoustics  can  be  used  to  evaluate 
the  phase  integral  producing  the  same  results  as 
direct  numerical  integration, 

-  show  that  the  two  methods  of  applying  ray  acoustics 
theory  (presented  in  Chapter  II)  produce  equivalent 
results  within  an  acceptable  tolerance, 

-  compare  speed  versus  accuracy  for  the  two  methods  for 
a  variety  of  test  cases,  and 

-  demonstrate  the  ability  to  identify  eigenrays  using 
either  of  the  ray  acoustics  methods. 


B.  DIRECT  INTEGRATION 

In  the  previous  chapter,  0y(y)  was  shown  as  the  solution 
to  the  phase  integral  in  the  depth  direction. 


©yCy)  =r''k  (C)dC 
^  Jyo  ^ 


(2.7) 


The  straight-forward  evaluation  technique  is  the  brute 
force  approach  using  direct  numerical  integration.  The  phase 
integral  will  be  evaluated  using  direct  numerical  integration 
for  a  simple  test  case.  The  results  will  then  be  compared  to 
phase  values  obtained  using  ray  acoustics  calculating  phase 
values  from  travel  times. 
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Conditions  for  the  test  case  are: 


-  source  depth  of  10  meters, 

-  receiver  depth  of  100  meters, 

-  linear  sound-speed  prpfile  with  a  single,  constant 
gradient  of  0.016  sec ^  and 

-  four  values  of  radial  spatial  frequency  (FR  values) 
will  be  evaluated. 

For  each  FR  value,  tabular  results  will  show  the 
corresponding  ray  launch  angle  The  relation  between 

the  FR  value  and  this  launch  angle  is  [Ref.  4]: 

/3(yo)  =  sin'^  [FR  *  c(yo)/f]  (3.1) 

where 

f  is  the  frequency  in  HZ,  and 

c(yo)  is  the  speed  of  sound  at  the  source  depth. 
Table  3.1  shows  the  results  produced  by  direct  numerical 
integration.  Tables  3.2  and  3.3  show  the  results  for  the 
Akima  cubic  spline/ordinary  differential  equation  solver 
(ODE  solver)  and  the  piecewise  linear  approximation  ray 
acoustics  techniques,  respectively. 

The  parameter  of  interest  is  THETAY  representing  the  phase 
change  in  the  depth  or  Y  direction.  Table  3.1  shows  these  0 
values  to  be  approximately  92,  86,  75,  and  56  radians  for  the 
four  FR  values  chosen.  Tables  3.2  and  3.3  contain  four  sets 
of  entries,  each  corresponding  to  an  FR  value.  These  entries 
begin  at  0  meters  range  and  end  at  the  range  corresponding  to 
100  meters  depth;  thus,  solving  the  phase  integral  for  the  0^ 
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TABLE  3.1 

PHASE  INTEGRAL  SOLUTION  BY  DIRECT  NUMERICAL  INTEGRATION 


INPUT  DATA  FOR  PHASE  INTEGRAL  EVALUATION 

F  •  2!,0.0  HZ  VL  ■  10.0  M  YU  ■  100.0  M  RATIO  •  NfR  •  S 

£ApS  •  0.000000  EREL  •  O.OOOICO 

YR£P  ■  0.0  M  CYREF  •  ISOO.O  M/SEC  C  ■  O.OUOOO  11/SEC) 

METHOD  •  2 


evaluation  by  IMSLlO  ROUTINE  DQOAG  WITH  IRULE  ■  2 


fricycles/m) 

BETAYO(DEG) 

RANGE  AT 

ANCLE  O'"  ARRIVAL 

TURNING  POINT 

TURNING  POINT 

THETAR(RAD) 

THETAV(RAO) 

' 

YU 

BETAYIOEG) 

DEPTHtH) 

RANCE(H) 

0.033216685* 

1 1.525 

16.3*009 

1 1.53* 

0.S7555E«0* 

0.65ie6E*0t 

0. 3e6oeE*oi 

0.*2291*79SiE*02 

0.01*5839718 

23.552 

39.25213 

23.57* 

0.16090E*0* 

0.21510E*0* 

0. 1*C2SE*D2 

0.e*S3860261E*02 

0.011863657* 

3*. 826 

*7.62812 

3*. 8*6 

0.62*82£*05 

0. 12522E»0* 

o.6:j:6e»o2 

0.75S7857  mE*02 

0.  i::i ;  M6*5 

52.069 

119.80722 

53.122 

0.23576E*0S 

0.70527E*05 

0. 10025E*OS 

0.5*573«8612e«02 

TABLE  3.2 

PHASE  INTEGRAL  SOLUTION  BY  ODE  SOLVER 

ERAVS  »  F  PRTALL  *  T 

YO  •  10.0  H  DEPTH  •  100.0  M  FREOC  •  250.0  H2  CHAX  •  1501.*  M/SCC 

FREOC/CMAM  •  0.14*<1  CYCLEj/H  ratio  •  FRHAX  •  0. 1***7  CYCLES/M  NFRS  •  5 

OLTFR  •  0.05!r96  CYCLES/H  RNCSTP  •  l.O  M  HRZRHO  •  150.0  M  VR  •  60.0  M 

MOTEi  FRMIN  •  OLTFR 


FR  •  Q.5S27660598W:*£-Ol  CYCLEO/M 

BETAO  •  n  ,5:6575211  709M  DEO 


rangE(H) 

CEPTH(M) 

TRVLT(CEC) 

BETAY(DEG) 

BTESr 

THETAT (RAD) 

THETAR(RAD) 

THETAY(RAD) 

0.0 

10.00 

0.000000 

11.5:5 

0. 155I8E-0J 

o.oooooooE*oo 

o.oooooooE*oo 

0.0000000E*00 

18.6 

100.00 

o.o*i:oo 

11.53* 

O.15J10E-OJ 

0.9tl3253£'02 

0.386086*E*01 

0.92291*8E*0: 

FR  •  o.**580«717*3<5:e-oi  cvclf:>/m 

BETAO  ■  rs.ssuoiBJfleou:  cco 


PANOC5M) 

DEPTHIM) 

TRVLTIOEC) 

BETayideg) 

btegt 

THETAT (RAD) 

THETARIRAD) 

TMETAY(RAD) 

0.0 

10.00 

o.ocoooo 

23.552 

O.2**3*E-03 

o.oooooooE*oo 

0.0000000£«00 

0.0000000E*00 

39.5 

100.00 

0.0*5620 

23.57* 

O.2**5*E-03 

0. 1027*1  1E>03 

0.1*62273E*02 

0.e*33860E*02 

FR  •  0. »9803657«651 78E-01  CVCLE5/M 

BETAO  •  S*.8:660An*66773  DEG 


PA»JGE  fM) 

DEP  TM ( M ) 

TflVLTfGEC ) 

BETAV(DEG) 

•  TEST 

THETAT(RAD) 

THETARIRAD) 

THETAY(RAO) 

0 . 0 

10.00 

0  -  000000 

T*  .826 

0. 31152E-03 

0  .OOOOOOOE‘00 

o.oooooooE*oo 

O.OOOOOOOE‘00 

*7.6 

100.00 

0.076132 

36 .86* 

0.591S2e-05 

0.  U  7  70:3E*03 

0.62S2376E*02 

0. 7537857E«02 

FR  *  0  .  1 23 1 7 !5:* ‘CO  CVCLES/M 

BETAO  •  53 . 06«:95  7582  1  3  15  DEG 

RA»K'EfMi  Dcr'TMfyj  TP,'LT($Er)  BETAVfrEG)  BTE5T  T»*ETATIRAD)  TMETAR(RADI  7HETAV(RAD) 

0.0  10. '^0  o.ccocco  52. O'-.’  o.5;r7iE-o3  0  .oooooooE*oo  0 .  occooopF  •  00  o.coocoonE«oo 

111.0  IOC. 00  0. 011001  55.172  O.S’ZTlE'OS  0  .  1 5 *8 2  * TE  •  0 5  0  .  I  00 Z 5 ?5E  •  0 5  0 . 5 *5  7 3 *8 E  •  0 2 
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TABLE  3.3 

PHASE  INTEGRAL  SOLUTION  BY  PIECEWISE  LINEAR  APPROXIMATION 


RAY  TRACING  PlElCEiRiC::  LIkCmR  APPRQKIMAT10^4 


KjATA  5 

Gt  1  )  = 

A 

0.  UOCCE-Cl 

l/GEC 

KE^'bER  OF  GRADIE^JTG  *  3 

G(  2)  »  O.UOOCE-Ol  l/GEC 

GC  3  )  = 

0. 16C0CE-CI 

1/CEC 

ERAVO  -  F  PRTALL  =  T 

VO  =  IC.C  H  DE^^TH  =  iCQ.O  M 

FPE'-T/Cf'A  •;  =  C.166A9  CYCLEC/M  RATIO  »  C.9'»09 
s  c.c::c^'.  cycles/m  rno^tp  ■  i.c  m 

^CT£:  FR‘^:n  =  L-TFR 


FR  «  C  .  3::  fpjo  I  ?r6E-Cl  C^CwEC/M 

EETAO  ■  U  .5:aS7SC  1  1  7053C  LEG 

RA‘,:E^'J  TRVwT(CLC)  EETAY(LEG)  BiEGT  ThETaT(RAD)  THETAR(RAD;  THETAVtRAO) 

io.jc  c.ogc.c30  o.i!jieE-o3  c.cccoocce^co  o.oocococe'OO  0 .  occcccoe*oc 

ia.-i  i:o.c:  c.coUl:  ii.sie  c.i::ieE-c3  c.9bi::siE*o:  o.iaAcs^eE*^!  o.9::9i68£*c: 


FR  =  c.-6:re?7;762<;s:E-:i  cvc.eo/m 

EE*a:  s  LEG 


ETEGT  ThETaTIRAD)  THETAR(RAD)  THETAVtRAD) 

o.:iie;t.E-c:  c.o:l:coce*co  o.ocooocoe*oo  o.coooocde*oo 

O.LfrbidE-c:  O.ic:7fell£*0:  0.16A2:73£*CL  0.86338‘iCE*02 


FR  ?  C  .  9?E3:  1  73E-C  1  CVClEG/M 

tlTAC  *  36.3LAAC9116AA7'6  L-3 

P-As:l.‘^)  TR.-TtEEE)  EETAviLEG)  B^EGT  TmETAT(PAD)  THETAR[RaL)  IhetayIRaD) 

1C. CO  C.CCCCCO  36.8CA  C.399S3E-C3  0  .  C  C  OCC' C  CE  •  0  0  0 . 0  C  C  0  0  0  C  E  •  C  0  0 . 0  C  C  0  00  CE  •  0  0 

t7.A  IOC. CO  C.C7‘-9::  I6.8b6  C.3?95:E-C3  C.n77CC3E‘C3  0  .  AC  30  3  7AE  ♦  OC  0 . 75  3  78  5  7E  ♦  CC 


»  : .  1 3 : ;  r  :  0  cyc^eo/h 

'AC  =  s:.c-^:ss;lu:i::o  leg 

'-.e:":  le“:h''*)  tr._7;ce:)  ee^avileg)  etegt  t.-e  r  at  )  thetapcr.^oj 

c.c  c.g':;:c:  si.cea  c.5::;ie-c:  c.ccc:cc:e*co  o.ccccooce*oo 

i:^.3  i::.c:  so.icc  c.55:7ie-c3  o.i5ta:6:E*c3  o.  ico:5:sE*c3 


TMET AV(PAC  ) 

0 . ccoocc:e*cc 

O.Sfc5;368E*OC 


LETAY(L'EG) 

C3.55: 

C  3 . 5  7  6 


FRECC  *  C5C.0  HZ  CM.^K  •  1501.6  M/OEC 

FR^x.^k  •  O.UbA;  CYCLEG/M  NURG  »  5 

HRZR».G  >  150.0  M  YR  a  AO.O  M 
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value  when  a  ray  travels  from  10  to  100  meters  depth.  The 
significance  of  these  three  tables  is  that  the  ray  acoustics 
methods  are  yielding  the  same  phase  values  0^  as  the  direct 
integration  method. 

The  conclusion  to  be  derived  from  this  is  that  ray 
acoustics  theory  can  be  used  to  accurately  evaluate  the  phase 
integral . 

C.  VARYING  SOUND-SPEED  GRADIENTS 

Both  ray  acoustics  methods  will  now  be  applied  to  three 
different  test  cases.  A  zero  gradient,  positive  gradient 
(+0.016  sec'^)  ,  and  a  negative  gradient  (-0.016  sec  linear 
sound-speed  profile  will  be  used.  These  simulations  will  show 
that  both  methods  produce  equivalent  results  within  an 
acceptable  tolerance  for  a  variety  of  media. 

The  conditions  common  to  each  simulation  run  are: 

-  source  depth  of  10  meters, 

-  horizontal  range  of  1  kilometer, 

-  ocean  depth  of  100  meters, 

-  speed  of  sound  at  the  surface  of  1500  meters  per 
second,  and 

-  five  values  of  radial  spatial  frequency  (FR) 
evaluated. 

The  results  are  presented  in  Tables  3.4  through  3.9.  The 
phase  integral  solution  is  the  phase  change  in  the  Y  (depth) 
direction  listed  as  the  MODULO  TWOPI  THETAY  value  in  radians. 
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TABLE  3.4 

ZERO  GRADIENT,  PIECEWISE  LINEAR  SOLUTION 


RAV  TRACING  USING  PICCENISC  LINEAR  APPROXIMATION 
CASE  lA 

NOATA  •  4  NUMBER  OP  GRADIENTS  •  S 

Cl  1)  ■  O.OOOOOE*00  t/SEC  Gl  2)  •  0.00000E*00  1/SCC 

G(  5)  •  O.OOOQOE*00  1/SEC 


VO  •  10. 0  M  depth  •  100.0  M 

FREOC/CMAX  ■  0.1444?  CYCLES/M  RATIO  • 

DLTPR  •  0.0335S0  CYCLES/M  RNCSTP  •  1.0  M 

NOTE:  PRMIN  >  DLTFR 


FREOC  •  2S0.0  HZ  CMAX  •  ISOO.O  H/SEC 

FRHAX  *  0.U44B  CVCLES/M  NFRS  •  5 

»«ZRNO  •  1000.0  M 


NFR 

BETAOIOEG) 

OEPTH(M) 

TRVLT(SEC) 

BETAYIDEG) 

THETATtRAD) 

THETAR(RAD) 

1 

1 1  .934 

90.51 

5.333447 

148.444 

e.9234SllE*04 

0.2094184E*0S 

2 

23.974 

98.44 

1 .444033 

194.424 

0.2418294E*04 

e.4188S71E«0S 

S 

34.844 

94.44 

1.111222 

143.134 

0.1749904C*04 

0.42829S7e*0S 

4 

93.  122 

39.79 

0.833417 

124.878 

0. 19O9120E*O4 

0.837474SE*05 

9 

09.190 

24.14 

0.444733 

89.190 

0.1047302E*04 

0. I047093E*04 

total  CPU  TIKE 

•  0  MIN 

.  22.45  SEC 

THETAVIRaD)  MODULO  TMOPI 
THETAV(RAO) 
0.502r0tSE*0«  0.S4499SI 

0.2199<!9E*04  0.50S7US 

0.U)7248E*04  9.124S041 

0.4714954E*05  0.2144804 

0.2094900E*0e  0.2094900 
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TABLE  3.6 

POSITIVE  GRADIENT,  PIECEWISE  LINEAR  SOLUTION 


RAY  TRACIMG  ySfW'J  PIECtWlSE  LINEAR  AfPROK |MAT IOW 
CASE  le 

NOATA  I  <  NUN8ER  OF  gradients  •  3 

G(  M  «  O.UOOOE-Ol  USEC  Gt  I)  ■  O.UOOOE-OJ  I/SEC 
G(  5)  •  O.liOOOE-01  1/SEC 

VO  •  10.0  M  DEPTH  •  100. 0  H  FREOC  •  250.0  HZ  CMAX  •  1501.4  M/SEC 

FPCOC/CHAX  •  0.|44<«  CVCtCS/M  RATIO  »  FRMAK  •  0.14447  CYCLES/M  HFRS  •  5 

DLTFR  •  0.033244  CVCLES/M  RNGSTP  •  1.0  M  HR^RWO  •  lOpO.O  N 

NOTE:  FRMIH  •  DLTFR 


NFP 

BETAO(DEG) 

DEPTMfM) 

TfiVLTtSEC) 

BETAY(DEG) 

TMETATIRAO) 

THETAR(RAD) 

1 

11.525 

97.80 

3.355452 

149.444 

0.5254499E‘06 

0.2091954C*05 

2 

25.552 

•4  .  *9 

1 .444924 

154.425 

0.241924)E*04 

0.4195908E*05 

5 

.  54  924 

55  39 

1.1111*2 

143.155 

0.1745454E^04 

0.4275843E*05 

4 

55.04* 

59 .73 

0  933379 

124.927 

0.  |509044E>04 

0.8S478I7E*03 

5 

97.343 

50.7  1 

0.447000 

97.974 

0. 1047721E‘04 

0.10459778*04 

total  CPU  time 

•  0  HIN 

.  22.05  SEC 

TMETAV(RAD) 

D.502729SE*06 
0.21M950E'04 
0.111 787DE*0C 
0.«7228«7E*0S 
0.l7<S79eE*01 


MODULO  TWOPl 
THETAVTRAD) 
0.  ?<<r549 
0 . 7S48**? 
5.7442109 
1  .0657792 
1.7457979 


TABLE  3.7 

POSITIVE  GRADIENT,  ODE  SOLVER  SOLUTION 


RAY  TRACING  USING  AKIMA  CUBIC  SPLINE  t  ODE  SOLVER 
CASE  29 

NDATA  •  4 

VO  •  10.0  M  DEPTH  •  100  0  H  FREOC  ■  250.0  HZ  CHAR  •  1501.4  M/SEC 

FREOC/CMAX  •  0.1444*  CYCLES/H  RATIO  •  0.**S9  FRMAX  •  0.14447  CVCLES/M  NFRS  •  5 

DLTFR  •  0.0352*4  CVCLCS'M  RNGSTP  •  |.0  M  HR2RNC  •  1000.0  M 

NOTE:  FRMtN  •  DLTFR 


NFR 

PETaOIDEG) 

DEPTH(M) 

TRVLTtSEC) 

BETAYfOEC) 

THCTATTRAD) 

THETAR(RAD) 

tmetayirad) 

MODULO  TWOPl 

TMETAV(RAD) 

1 

n  .525 

97.80 

5 .555452 

149.444 

0.5r54499C*04 

e.2091954E*0S 

0.50r7295E*04 

0. 74424*5 

2 

25.552 

94. *9 

1  . 44«9:4 

154.425 

o.r*ii24i€*o4 

0.*lt5*08E*05 

0.ri**'950F»?< 

0  .  r 3**005 

5 

34.824 

55.58 

1,1111*2 

143  155 

0.  |745454E*C4 

0  4275945E‘05 

0. Ill ?870E*04 

5.7442125 

4 

55.049 

59 . 7  5 

0.«5:5TS 

124.927 

e.l509044E*04 

O.B547817E*OS 

o 

CA 

O 

1.04577*5 

5 

97.343 

50.71 

0 .447000 

97.974 

0. I047721E*04 

0.1045977E‘0* 

0.1743799E*01 

I  .74J  79  79 

total  CPU  TIME  •  14  him  .  r<.50  SEC 
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TABLE  3.8 

NEGATIVE  GRADIENT,  PIECEWISE  LINEAR  SOLUTION 


RAY  TRACING  USING  PIECEWISE  LINEAR  APPROXIMATION 
CACE  1C 

NCATA  »  A  NUMBER  OF  GRADIENTS  •  5 

C(  I)  •  -0.I<000E-01  i/SEC  GC  Z)  *  -O.IAOOGC-01  I/SEC 
G(  S)  «  -0.  UOOOE-Ol  l/SEC 

YO  •  10. 0  M  depth  »  100. 0  M  FREOC  •  2S0.0  HZ  CHAX  ■  1500.0  M/SEC 

FRECC/CMAX  •  O.U4t7  CVCLES/M  RATIO  •  0.44t4  FRHAX  •  0. 14445  CYCLES/M  NFRS  •  5 

Ol^FR  •  O.OJ!350  CYCLES/M  RNG5TP  •  1.0  M  MRZRNC  •  1000.0  M 

NOTE;  FRMIN  ■  DLTFR 


NFR 

BETAOTOEG3 

OEPTHIM) 

TRvlTISEC) 

BETAVIDEC) 

THETATIRAD) 

THE  TAR (RAD) 

TMETAY(RAO) 

MODULO  TWOPI 

THEIAy(RAO) 

1 

11  .535 

87 . 78 

5 .  !372<0 

148. A75 

0.5242124E*04 

0.2094ie4E*05 

0.5052706E*04 

4.1574514 

2 

23.57! 

94  .  48 

1 .468422 

154.450 

0.242I045E*04 

0.4188!71E*05 

0.2202228E«04 

S. 11  34455 

3 

34.861 

55 . 32 

1  .  1  12^:8 

143 . 140 

0.1  747414E*04 

0.4282557E»0S 

0.  1U9159E*04 

0.7515715 

« 

53  .  1  U 

38.4  3 

0  .eSA»AS 

124 , 909 

0.  t3l0586E«04 

0.8S74743E*0S 

0.4729117E*0S 

1.4727921 

5 

88.835 

35  .  *7 

0.6470A6 

88.224 

0. 1047793E‘04 

0 . 104709SE*04 

0.4999412E*00 

0.4994412 

total  CPU  time  •  0  MIN  .  :z.4l  SEC 


TABLE  3.9 

NEGATIVE  GRADIENT,  ODE  SOLVER  SOLUTION 


RAY  tracing  using  AKIMa  CUBIC  SPLINE  I  ODE  SOLVER 
Case  :c 

NDATA  «  « 

FREOC  •  250.0  MZ  CMAX  •  1500.0  M/SEC 

FRMAX  •  0.14445  CYCLES/M  NFRS  >  5 
HRZRMG  •  1000.0  M 


▼0  »  lO.O  M  DEPTH  •  100.0  H 

FOEOC/CMAX  •  0.14447  CvClES/M  RATIO  •  0.*44» 
OLTFfl  ■  0.055330  CYClES/M  RNGSTP  •  1.0  M 

NOTE:  FRMIH  •  DLTFR 


NPR 

BETAO(OEG) 

OEPTMIM) 

TRVLT(SEC) 

BETAV(0€C) 

THETAT(HAD) 

TM£TAR(RAD) 

1 

1  1  .535 

87.77 

3 .337240 

168.475 

0.5242124C«04 

0.2094I84E*05 

2 

23  5  73 

^6  .*8 

1  . 648622 

156.450 

0.2421045C*04 

0.4te8S7lE*0S 

3 

34.861 

55 . 32 

1  .  1  12438 

143.160 

0.  I747414E«04 

0.4282557E*05 

A 

53.114 

38.65 

0.834345 

124.909 

0.  I3)0584E«04 

0.857S76XE»0! 

5 

88  8!5 

35.67 

0 . 46  7946 

88 .224 

0.  I047793E»04 

0. 104709!E*04 

thetav(rad)  modulo  twopi 

TMETAY(RAD) 
0.503Z704E*0<  4.1576454 
0.2202228E»0<  S. I15<4e3 
0.mfl5«E*0<  0. 75157:4 
0.A72»W7E«05  1.4727421 
0.444«61:E»00  0.444»41: 


total  CPU  time  •  14  MIN  .  21.24  SEC 
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Tables  3.4  and  3.5  show  the  results  for  the  zero  gradient 
case.  A  comparison  between  two  tables  shows  that  the  phase 
values  (MODULO  TWOPI  THETAY)  agree  perfectly.  Additionally, 
all  other  calculated  values  in  the  tables  agree.  The  only 
difference  is  the  TOTAL  CPU  TIME.  Both  methods  are  using  a 
conservative  range  step  (RNGSTP)  of  1  meter,  i.e.,  the 
incremental  increase  in  range  is  1  meter  in  the  propagation 
calculations.  For  these  conditions,  the  noE  solver  is  taking 
over  30  times  longer  to  run  the  simulation. 

Tables  3.6  and  3.7  show  the  results  for  the  positive 
gradient  case.  A  comparison  between  these  two  tables,  as 
well  as  the  negative  gradient  results  shown  in  Tables  3.8  and 
3.9,  shows  acceptable  agreement  between  the  two  very 
different  calculation  methods.  The  only  discrepancies  in 
calculated  values  occurs  in  the  modulo  27r  0^  phase  values 
needed  for  solving  the  phase  integral.  The  phase  calculation 
is  most  challenging  because  of  the  nature  of  phase  (being  a 
modulo  27r  function)  .  The  discrepancies  for  both  the  positive 
and  negative  gradients  are  on  the  order  of  hundred-thousandths 
of  a  radian  or  less.  The  CPU  times  continue  to  follow  the 
pattern  seen  earlier  that  the  ODE  solver  requires  over  3  0 
times  longer  to  complete  a  simulation  run. 

The  following  ray  trace  plots.  Figures  3.1  through  3.6, 
correspond  to  Tables  3.4  through  3.9.  The  plots  assist 
interpretation  of  the  data,  but  as  expected  from  the  tabular 


20 


QDOQoy  r  f  HNJ 


hRZHNG:  1.0  KM  TR:  ‘10.0 


RflY  TRRCE:  PIECCtjlSE  lUJERR  “^PPROXIMRIION 


!0.0  M  HRZRNG: 


PftY  TfflCC:  FlKira  CUBIC  SPLINE  I  ODE  SOLVER 


HRZRNG: 


ROY  TRFICC:  P'CCCl^ISc:  LINCRR  RPPROX  IWT  1  ON 


HRZRNG:  1.0  KM 


results,  no  variations  in  depth  or  angle  of  arrival  can  be 
detected. 

D.  CHALLENGING  TEST  CASES 

A  smooth  sound-speed  profile  in  the  form  of  a  half 
sinusoid  describes  the  medium  for  the  next  set  of 
simulations.  To  this  point,  the  agreement  in  accuracy  between 
the  two  methods  has  been  excellent;  however,  the  shape  of  a 
sine  curve  is  difficult  to  approximate  with  piecewise  linear 
segments.  While  this  will  be  a  challenging  test,  smoothly 
curved  profiles  do  occur  naturally.  The  following  set  of 
values  were  used: 

-  1500  m/sec  is  the  sound  speed  at  the  ocean  surface  and 
at  the  2000  meter  bottom,  and 

-  a  minimum  sound  speed  of  1475  m/sec  occurs  at  the  1000 
meter  depth. 

Three  sets  of  simulations  are  presented  using  5,  11,  and 
17  sound-speed  versus  depth  data  pairs.  The  simulation 
results  will  be  examined  to  quantify  the  number  of  ocean 
medium  sound-speed  samples  required  by  each  method  to  converge 
to  a  solution.  Also,  the  agreement  in  results  for  the  two 
methods  and  the  CPU  times  will  be  evaluated. 

As  the  number  of  equally  spaced  data  pairs  varies,  the 
medium  and  source  conditions  will  be: 

-  source  depth  of  1000  meters, 

-  horizontal  range  of  3  kilometers,  and 

five  values  of  radial  spatial  frequency  (FR) 
evaluated . 
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Tables  3.10,  3.11,  and  3.12  show  the  results  for  the  5, 
11,  and  17  data  pair  test  cases  for  the  piecewise  linear 
method.  In  comparing  the  first  two  of  these  tables,  large 
discrepancies  are  seen  in  the  final  depths  for  each  of  the 
five  radial  frequencies  evaluated — errors  as  large  as  seven 
meters.  Likewise,  travel  time  values  differ  by  hundredths  of 
a  second  meaning  that  modulo  2tt  phase  values  show  no 
correlation  between  test  cases.  Unfortunately,  the  variations 
between  Tables  3.11  and  3.12  where  more  data  points  are  used, 
show  the  piecewise  linear  method  is  not  converging  to  a 
solution . 

Figures  3.7  through  3.9  correspond  to  Tables  3.10  through 
3.12.  While  numerical  phase  values  show  gross  disagreement, 
the  inability  to  detect  differences  in  the  graphical 
representations  illustrates  the  very  sensitive  nature  of  the 
phase  calculations. 

Tables  3.13  through  3.15  show  the  results  when  the  ODE 
solver  runs  the  simulations  for  5,  11,  and  17  data  pairs.  The 
5  data  pair  case  of  Table  3.13  again  shows  significant 
disagreement  with  the  Table  3.14  results  using  11  data  pairs. 
Final  depths  show  differences  of  as  much  as  a  meter  while 
travel  times  vary  by  milliseconds.  As  with  the  piecewise 
linear  method,  the  modulo  27r  phase  values  have  not  yet 
converged  to  a  solution.  Variations  of  up  to  2  radians  are 
seen . 
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TABLE  3.10 

FIVE  DATA  PAIRS,  PIECEWISE  LINEAR  SOLUTION 


RAV  TRACING  USING  RieCEWlSE  LINEAR  APPROK IMAT ION 
CASE  ISOS 

NOATA  •  s  NUMBER  OF  GRADIENTS  •  < 

G(  I)  •  -0.5S3S5E-01  l/SEC  C(  ■  -0.IR4ASE-01  I/SEC 
C(  5)  •  0.IAA<SE-01  l/SEC  G(  AJ  •  0.5SSSSE-01  l/SEC 

ERA/S  •  F  PRTALL  •  F 

VO  •  1000. 0  M  depth  •  2000.0  M  FREOC  •  2S0.0  HZ  CMAK  •  1500.0  M/SEC 

FREOC/CMAK  •  0.14««7  CVCLES'H  RATIO  •  0 - FRMAX  •  0.144*5  CVCLES/H  MFRS  •  5 


DLTFfl  •  0.033330  CvcLES/M 

NOTE;  FRMIN  •  DLTFR 

RNGSTP  • 

1 .0  H 

HR2RHG  • 

5000. 0  K 

YR  • 

AO.O  M 

K^R 

GETAOfDEG) 

DEPTH(M) 

TRVLTfSEC) 

BETAYfOEG) 

THETAT(RAD) 

THETAR(RAO) 

TMETAYfRAD) 

MODULO  THOP; 

THETAY(RAD) 

1 

11.3A1 

lAA . 7A 

10.207184 

148. 50A 

O.UOSSA1E*05 

0.4282557E*0S 

0.15A05I5E«05 

5.047A4AS 

2 

23.140 

A1  .81 

5.  1  033  1  1 

154. AA4 

0.8014242E*0A 

8.1254511E*0A 

0.475975IE*OA 

5  •  5244822 

5 

34.153 

104A .44 

3.A0:a34 

34. 180 

O.S3AASS8E*04 

0. 1S8A747E*0A 

0.5A54771E*0A 

*‘.0141220 

A 

51  .B46 

482  .A7 

2  55A7A7 

127.402 

O.A012«87E«OA 

0.25ISC23E*0A 

0. |A4444AE*0A 

4.5458515 

5 

74. A** 

1510.35 

2,057830 

81 .2A4 

O.S252A53£*OA 

0. 5U1278E*0A 

0.4115AOAE*02 

5.184AA04 

total  CPU  TIME  •  0  MIN  .  2A.A5  SEC 


TABLE  3.11 

ELEVEN  DATA  PAIRS,  PIECEWISE  LINEAR  SOLUTION 


RAV  TRACING  USING  PIECEnJSE  LINEAR  APPROK I  MAT  ION 
CASE  ISU 

NOATA  •  11  NUMBER  OF  GRADIENTS  •  10 

C(  1)  •  -0.3t427£-0l  1 'SEC  G(  O  •  -0. SABAAE-OI  l/SEC 

G(  S)  •  -0.27«SAe-0l  ’/SEC  G(  A|  «  -0.I7755E-01  |/SEC 

G(  5)  •  -O.4M7«£-02  l/SEC  Of  4)  •  0-41ir«E-0r  l/SEC 

G(  7)  •  0.17755E-01  l/SEC  Of  B)  •  0.2;45AE-0I  l/SEC 

Of  9)  •  O.JA8A4E-01  l/SEC  GflO)  •  0.JB427E-01  l/SEC 

ERAYS  •  F  PRTALL  •  F 

VO  ■  1000,0  M  depth  •  2000.8  M  FREOC  •  250.0  M2  CMAX  ■  1500.0  M/SEC 

FREOC/CMAX  •  0.1  444  7  C^CLES/M  RATJO  •  O.V***  FRMAX  •  0.14445  C'<CLES/H  MFRS  •  5 

DLTFfl  •  0.CS5J50  CvCLES'm  RNGSTP  •  l.O  M  MRZRMG  •  5000. 0  M  VR  •  *0.0  M 

NOTE;  FRHJN  ■  DLTFfl 


NFR 

BETAOfDEG  J 

DE®tn(mi 

TRVLTISEC) 

iETAY(DEC) 

THETATIRAD) 

THETaR (RAD) 

TMETAY(RAD) 

MODULO  TH(5P1 

TMETAY(RAD) 

1 

1 1  T<  1 

157.50 

10. 

214744 

148.504 

0. 14Q4844E*05 

0.i:8:557E»05 

0.\54:020E*05 

1 .2428*44 

2 

23.140 

37. *7 

5  . 

108375 

154.4*4 

8.80257*rE»0< 

0  i:545nE*04 

0.4747231E*0< 

0.2400057 

S 

34.  153 

104  7  5* 

3  . 

405626 

54.  165 

0  5J<*545E*0< 

0- 188C747E*D4 

0.I64*778E*04 

2 . 7450820 

A 

51  .868 

47«.33 

2  . 

5575*4 

127.445 

0.4M7070E*04 

0.2S130:3E‘04 

0. 150«048E*0* 

2.36436«1 

5 

74 . 6«A 

1520.40 

040166 

81 . 354 

0. 3234104E»0« 

0. 5141Z78E*0* 

0.4*e:580E*02 

0.5780145 

total  CPU  time 

•  0  MIN 

.  25 

.34  SEC 
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TABLE  3.12 

SEVENTEEN  DATA  PAIRS,  PIECEWISE  LINEAR  SOLUTION 


RAY  TRACING  USING  PIECEWISE  LINEAR  APPROXIHAT ION 
CASE  isi; 


NDATA  •  17 

Gin-  -0.5«0I8E-0I  I/SEC 
G(  SI  ■  -0.S4577E-01  I/SEC 
G(  51  •  -0.2<i73E-01  I/SEC 
GI  71  •  -0. 11581E-01  1/SEC 
G(  SI  •  0.38«rSE-02  1/SEC 
Gtlll  •  0.18<82E-01  1/SEC 
OllSl  •  0.SOSO7E-01  1/SEC 
OI15I  .  0.3751SE-01  1/SEC 


NUMBER  OP  GRADIENTS  •  I( 

G(  21  •  -0.S751SE-01  1 /SEC 
G(  «1  ■  -0.S0SO7E-O1  1/SEC 
G(  «l  •  -0.ie<82E-ei  1/SEC 
G(  81  •  -0.!e«2«E-02  1/SEC 
GMOl  •  0.11S81E-01  1/SEC 
G(121  •  0.2A873E-0I  1/SEC 
0(1«1  •  0.S<B77E-01  1/SEC 
0(141  •  0.SS018E-01  I/SEC 


ERAVS  •  F  PRTALL  •  F 

YO  ■  lOOO.O  M  DEPTH  .  2000.0  M  FREOC  •  250.0  H2  CMAX  •  1500,0  M/SEC 

FREOC/CMAX  ■  0.14447  CYCLES/M  RATIO  •  O.SSSS  FRMAX  •  O.U(45  CYCLES/M  NFRS  •  5 

OLTFR  •  O.OSSSSO  CYCLES/M  RNGSTP  •  1.0  M  HRZRNO  •  5000.0  M  YR  •  40.0  M 

NOTE!  FRMIN  .  DLTFR 


NFR 

betaoidegi 

OEPTHIMl 

TRVLTISECl 

BETAY(DEG) 

TMETATIRADl 

THETARIRAD) 

THETAVIRADl 

MODULO  THOPl 

THETAVIRADl 

1 

ll.SAl 

154.47 

10.217844 

148.504 

0.14050I»E»05 

0.4282557E>05 

0.1542U5E*05 

2.5957577 

2 

25.140 

57.55 

5.108420 

154. 44S 

0.8024402E*04 

0. 12545I1E«04 

0.47480SIE*04 

1 .1002258 

S 

54.155 

1047.73 

5.405SS4 

54.140 

0.5350I24E«04 

0.1884747E*04 

0.544535»E*04 

5.5242727 

4 

51.848 

4  78.S8 

2.557444 

127. S75 

0.«0I7558E*0« 

0.25I5025E*0« 

0.15045I5E*04 

2.8558755 

5 

7S.4S4 

15.':  .42 

2.040452 

81.522 

0.5254550E*04 

0.5141278E*04 

0.»527150E*02 

1.0257247 

TOTAL  CPU  TIME  •  0  MIN  ,  25 . S<  SEC 
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CASE:  1S05 

YO:  1000. J  M  HRZRNG:  3.0  KM 

RRYS:5 


’iCOO.Q  HRZRNG:  3.0  KM 


TABLE  3.13 

FIVE  DATA  PAIRS,  ODE  SOLVER 


RAY  tRACINO  USING  AK|MA  CUBIC  SPLINE  •  ODE  SO(.VCR 
CASE  2S0S 

NDATA  •  S 

ERAYS  •  E 

VC  •  lOOC.O  M 

^BEOC/CHAR  •  0.141*7  CVCLES/M 
DLTPR  ■  O.OJJJSO  CVCLES/H 
NOIC;  FRHIH  •  DLTPR 


MTR 

BETAOIOEGI 

OCPTHtHJ 

TRVLTISECl 

fETAvtOEGI 

TMETATIRAOI 

THETARIRADI 

THCTAV(RAD) 

HODULO  THOR] 

THETAVtRADl 

J 

JI.34I 

154.05 

10-220257 

140  50* 

0. 1405344E*05 

e.42S2557£<0S 

e.l542544e*05 

0.4444177 

2 

2!.  140 

34  54 

5  .  I0*8:s 

154  452 

e.80r4441E*0« 

0.12545lie*04 

0.i744480E*04 

2.4B4I844 

3 

34  .  153 

1048  44 

3.404808 

34.157 

e.5S5i4orE*o« 

0.  IB84747E*04 

e.S4444S5E*04 

4  .5447015 

« 

51  .048 

478.30 

2.558144 

127.473 

O.40t8402e*04 

0.2: is025e*04 

0. 1505574e*0< 

1.4470735 

5 

74  444 

1521 .82 

2 .0*04«8 

81.304 

0.S2S4425E*04 

0.3141278E*04 

0. 4SS4448E*02 

1.0440 400 

RRTALL  •  F 

DEMH  •  2000.0  M  PREOC  •  2S0.0  K?  CHAM  •  tSOO.O  N/$EC 

RATIO  •  0.9*«*  fRHAX  •  CYCLES/H  NFRS  •  $ 

RMGSTP  •  1-0  M  HRZRNC  •  S000.9  M  VR  •  40.0  M 


TABLE  3.14 

ELEVEN  DATA  PAIRS,  ODE  SOLVER 

RAY  TRACINO  USING  ARIH4  CUBIC  SPLINE  I  OOC  SOLVER 
CASE  2Sn 

MOATA  •  11 

fRAVS  ■  f  RRTALL  •  E 

VO  •  1000. 0  H  CEPIM  •  2009.0  M  FRCOC  •  2S0.0  HZ  CHAM  •  ISOO.O  M/SEC 

RBCOC/CMA*  •  0.1  444  7  CyCLES-'M  RATIO  •  0.t«44  ERKAK  •  0.14441  CVCLES/M  NTRS  •  S 

OLTfR  •  0.0S5J1C  CVCLES/H  RnGSTP  •  1.0  M  HRZRNC  •  3000.0  M  VR  •  40.0  M 

NOTE«  RRMIM  •  OLTT-R 


MFR 

BCTAQIOCO) 

DtPTMIMl 

TttVLTtsccJ 

ttTAVIDECl 

TTCTATfRADl 

THETAR(RAD) 

TKTAVIRAD) 

MODULO  TNOPI 

THCTAVIRAD) 

1 

11.341 

154.04 

10.218425 

148.504 

0.|*05!I«E»05 

0.ir8Z557E*eS 

0.15423121*05 

4.  I8Z8«BZ 

2 

23.  140 

57.23 

5.  108*4; 

154  444 

0.80Z5144C*04 

0.  l2545nE*04 

e.4748i8SE*0« 

1  .44Z04Z5 

3 

54.153 

1047.47 

5.404254 

54.157 

0  5350551E»04 

0.18847I7E*04 

0.3445744E*0« 

3.7:88771 

4 

51  .848 

478.73 

2.557848 

127.470 

0.4CI7858E‘04 

e.Z5l5025E‘04 

0.1504835E*04 

3.1540545 

5 

74.444 

1522  ST 

2.040455 

tl  .312 

0.SZ34B34E*04 

0.S14I278E*0« 

0.4555554E*02 

1.3078080 

TOTAL  CPU  time  •  <0  HIM  .  M.SJ  SEC 

TABLE  3.15 

SEVENTEEN  DATA  PAIRS,  ODE  SOLVER 

RAV  TRACING  USING  AYIMA  CUBIC  SPLINE  I  OPE  SOLVER 
CASE  2S17 

MOATA  •  17 

ERAvS  •  W  PSTALL  •  P 

VO  •  10C0.O  H  PEPTM  •  ZOOO.O  M 

ERECC^’CHA*  •  0.14S4?  CvCLES'H  BATJO  •  0.S*T* 

DLfPR  •  0  0!!St0  CYCLES'**  HwGgtp  •  1.0  M 

NOTE;  FRMjw  ,  dltfr 


EREOC  •  2S0.0  HZ  CHAM  *  1500.0  H/SEC 

FB**AM  •  0.14445  CVCLES/M  NPRS  •  5 

MRZDMG  •  3000.0  M  VR  •  40.0  H 


MOPLLO  TMOPl 
TMETAY(RAD) 
4.1I54B74 
1.4570417 
3.  706B54t 
3.IJ5S514 
I  ?T545Bf 

TO’AL  CPU  time  .  *;  NJW  ,  SEC 

14 


WPR  BETAOIPESI  0€''TMfM)  T9VLT(SECT  BE’AvfPEGT  THETATfRADT  TMCTARtHAO)  T>«TAV(RAD) 

1  n  134. IJ  15  ZlfSRO  14B.55*  0 . ) »09 1  3  IE *05  0 . 4ZBZ55  7E  •  C 5  0 . 15<2 505E • 05 

2  2  5  1  40  5  7-25  5  1  08  475  154  . 0 .0025 1 54E  •  04  0.IZ545nE*04  0 . 474t  4C8E  *  0« 

5  54  15  5  IC47.44  5,CC«:S*  S4.l'7  0 . 5 55050 7E ‘O*  0  .  IBB4 74 7E • 04  0 . 54  4  5 740E • 04 

4  51.8*8  *78  7<  Z.55'8I4  Ir7-»7f  0 . <0 t ;BS«E *04  0 . Z5 1 502 3E • 04  0 . 1 5 0«8 1 7E • 04 

5  74.<4C  1  522  Ji  2-045425  #l.!l  7  0 . 32 5482  2E • 04  0 . 3 J 4 J 2  78E  •  04  0 . 4554 5 ?4E *02 


A  comparison  of  Tables  3.14  and  3.15  shows  that  phase 
calculation  results  converged  for  each  of  the  five  FR  values. 
The  increase  in  data  pairs  to  17  from  11  has  changed  the 
final  depth  values  by  only  a  few  hundredths  of  a  meter. 
Likewise,  the  very  sensitive  modulo  2"^  phase  values  show  a 
change  of  only  hundredths  of  a  radian  or  less.  The  ray  trace 
plots  for  Tables  3.13  through  3.15  are  provided  in  Figures 
3.10  through  3.12. 

The  CPU  times  for  these  three  tables  show  the  ODE  solver 
to  be  costly,  but  insensitive  to  increases  in  the  number  of 
data  pairs.  That  is,  approximately  the  same  amount  of  CPU 
time  is  reguired  regardless  of  the  number  of  data  pairs  used. 
The  CPU  time  for  the  piecewise  linear  method  also  proved 
fairly  insensitive  to  the  number  of  data  pairs  used  as  seen  in 
Tables  3.10  through  3.12. 

The  significant  findings  of  this  section  are: 

-  the  ODE  solver  can  perform  accurate  phase 
calculations  with  only  17  data  pairs  sampled  from  a 
2000  meter  deep  ocean  having  a  curved,  sinusoidal 
sound-speed  profile, 

-  the  piecewise  linear  approach  does  not  converge  to  a 
solution  using  the  17  data  pairs,  and 

-  the  differential  equation  solver  is  very  costly  to  use 

in  terms  of  CPU  time  versus  the  piecewise  linear  approach. 


35 


RPY  TRACE:  '^Kim  cuaic  SPLINE  I  ODE  SOLVER 


CHSE:  2S05 

TO:  1000.0  M  HRZRNG:  3.0  KM 
RAYS: 5  _ 


PKlMfl  C'JSiC  5^L1Nl  I  ODl  solver 


HRZRNG;  3.0  KM 


PPY  TI^RCC:  FIKlMFI  CUBIC  5*^1  !Nl  i 


HR2RNG:  3.0  KM 


E.  SPECIAL  INVESTIGATIONS 


Each  of  the  two  phase  calculation  methods  has  shown  one 
characteristic  weakness.  The  piecewise  linear  method  did  not 
converge  to  a  solution  in  the  sinusoidal  sound-speed  profile 
case  when  17  data  pairs  or  samples  were  provided.  In 
contrast,  the  ODE  solver  does  converge  to  a  solution,  but  is 
costly  in  terms  of  CPU  time.  This  section  will  attempt  to 
tailor  the  settings  for  each  method  to  compensate  for  its 
weaknesses . 

The  simulation  results  in  Tables  3.16  through  3.18  were 
produced  by  the  piecewise  linear  method  using  29,  55,  and  65 
data  pairs,  respectively.  All  medium  and  sound-speed  profile 
characteristics  remain  the  same  as  in  the  previous  section. 
As  the  number  of  data  pairs  increases,  the  6^  values  are 
converging  to  the  Table  3.15  solutions  obtained  with  the  ODE 
solver.  Unfortunately,  the  phase  values  continue  to  show 
unacceptable  errors  of  tenths  and  hundredths  of  a  radian  for 
the  best  case  of  65  data  pairs. 

Table  3.19  shows  the  results  for  the  ODE  solver  when  the 
range  step  (RNGSTP)  parameter  is  increased  to  five  meters.  As 
explained  earlier,  the  range  step  is  the  differential 
equation  system  independent  variable.  All  calculations  prior 
to  Table  3.19  used  a  conservative  range  step  of  one  meter  for 
both  methods.  Table  3.19  compares  favorably  with  Table  3.15 
achieving  a  balance  between  CPU  time  and  accuracy.  While  the 
degradation  in  accuracy  is  only  microseconds  for  travel  time 
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TABLE  3.16 

TWENTY-NINE  DATA  PAIRS,  PIECEWISE  LINEAR  SOLUTION 


RAV  TRACINO  USING  PIECEWISE  LINEAR  APPROXIMATION 
CASE  1S2T 


NDATA 

• 

25 

NUMBER 

OF  GRADIENTS  ■ 

28 

G(  1  I 

• 

-0.3yi8SE-0l 

I /SEC 

01  2) 

-0.58<55E-OI 

1/SEC 

01  5) 

• 

-0.577I5E-01 

1/SEC 

0(  <) 

-O.S«2*2E-OI 

1/SEC 

G(  5) 

• 

-0.5<U52E-01 

I /SEC 

0(  i) 

-0.32010E-0I 

1/SEC 

Gl  7) 

■ 

-0.252i«E-01 

1/SEC 

0(  8) 

-0.2AIS4E-01 

1/SEC 

01  ») 

• 

-0.227I2E-01 

1/SEC 

GUO) 

-0. 1858*E-01 

1/SEC 

GUI) 

• 

-0.15020E-01 

1/SEC 

GU2) 

-0. 108<«E-01 

1/SEC 

GUS) 

■ 

-0.»57«5E-02 

1/SEC 

G(  14  ) 

-0.22007E-02 

1/SEC 

0(15) 

■ 

0.22007E-02 

1/SEC 

out) 

O.IS745E-02 

1/SEC 

GU7) 

■ 

0. lOOtiE'Ol 

1/SEC 

0(18) 

0.15020E-01 

1/SEC 

GUS) 

■ 

0. 1858<E-01 

1/SEC 

G(20) 

0.22712E-0I 

1/SEC 

0(21) 

• 

0.2tl5AE'01 

1/SEC 

0(22) 

0.252««E-01 

1/SEC 

0(25) 

m 

0.520t0E>01 

1/SEC 

0(24) 

0.54352E-01 

1/SEC 

0(25) 

a 

0.S(2«2E'01 

1/SEC 

G(2i) 

0.57715E-01 

I /SEC 

0(27) 

a 

9.58t5SE-01 

1/SEC 

0(28) 

0.35188E-0I 

I/SEC 

ERAYS 

a 

F 

PR TALL  •  F 

YO  •  1000.0  M  DEPTH  .  2000.0  M  PREOC  •  250.0  KT  CHAR  .  1500.0  M/8EC 

FREOC/CMAX  •  O.Uttr  CYCLES/M  RATIO  •  O.TPTt  FRMAK  •  O.UtiS  CYCLES/M  MFRS  •  5 

DLTFR  •  O.OSJJJO  CYaES/M  RNCSTP  •  I.O  M  HRZRNO  •  5000.0  N  V«  •  <0.0  H 

NOTEl  FRMIN  •  DLTFR 


NFR 

7 

BETAO(DEG) 

DEPTH(M) 

TRVLT(SEC) 

BETAY(DEO) 

THETAT(RAD) 

THETARIRAD) 

THETAY(RAD) 

MODULO  THOFI 

THETAY(RAD) 

1 

11.541 

154.51 

10.218542 

148.504 

0.1405055E*05 

0.4282557E>D5 

D.1542248E*05 

5.7415052 

2 

25.140 

57.54 

5. 1088S7 

154.445 

0.8024575E*04 

0.1254511E*04 

0.4748442E*04 

1.4710844 

5 

54.155 

1047.88 

5.404158 

54.157 

0.55505e0E*04 

0.18B4747E*04 

0.544S415E>04 

5.5775585 

4 

51.848 

478.82 

2.557772 

127.574 

0.4017755E*04 

0.25I5025E*04 

0.1504714E<04 

5.0548745 

5 

7».4»4 

1522.11 

2.040547 

81.520 

0.5254750E*04 

0.514I27BE«04 

0.5545157E*02 

1.2041554 

TOTAL  CPU  TIME  ■  0  HIM  ,  2i.<«  SEC 
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TABLE  3.17 

FIFTY-FIVE  DATA  PAIRS,  PIECEWISE  LINEAR  SOLUTION 


RAY  TRACINO  irSlNG  PIECEMISE  LINEAR  APPROXIMATION 
CASE  1SS5 


NDATA 

5$ 

NUMBER  OF  gradients  • 

54 

G(  n 

-0.33248E-01 

1/SEC 

G(  2) 

-0.33115E-01 

1/SEC 

G(  5) 

-0.38850E-01 

1/SEC 

G(  4) 

-0.384S3E-01 

1/SEC 

G(  S) 

-0.3732tE-0I 

1/SEC 

G(  t) 

-0.37271E-01 

1/SEC 

G(  7) 

-0.34430E-01 

1/SEC 

0(  8) 

-0. 3$5e4E-01 

1/SEC 

0(  ») 

-0.34$4IE-ei 

1/SEC 

0(10) 

-0.334I3E-0I 

l/SEC 

GUI) 

-0.321i3E-0I 

1/SEC 

0(12) 

-0.30800E-01 

1/SEC 

G(I3> 

-0.2333IE-01 

I/SEC 

0(14) 

-0.27744E-01 

1/SEC 

G(IS) 

-0.24103E-01 

1/SEC 

0(14) 

-0.24353E-01 

1/SEC 

0(17) 

-0.22521E-01 

1/SEC 

G(18) 

-0.20413E-01 

l/SEC 

G(l«> 

-0. 1843$E-0l 

1/SEC 

0(20) 

-0. 14$34E-01 

1/SEC 

G(21 1 

-0. 14437E-01 

1/SEC 

G(22) 

-0. 12350E-01 

1/SEC 

0(23) 

-0.10I42E-0I 

1/SEC 

0(24) 

-0. 73333E-02 

l/SEC 

0(2$) 

-0.54307E-02 

1/SEC 

G(24) 

-0.34221E-02 

1/SEC 

0(27) 

-0.11«20E-a2 

1/SEC 

0(28) 

0.11420E-02 

1/SEC 

0(23) 

a.34221E-02 

1/SEC 

0(30) 

0.54307E-02 

l/SEC 

G(31) 

0.73333C-02 

1/SEC 

0(32) 

0. 10142E-01 

1/SEC 

0(33) 

0. 12S50E-dl 

1/sec 

0(34) 

0. I4437E-0I 

l/SEC 

0(5$) 

0.  U$3«E-01 

1/SEC 

0(34) 

0.I843$C-01 

l/SEC 

0(37) 

6.20(13E-01 

1/SEC 

0(38  ) 

0.22$21E-01 

l/SEC 

0(33) 

0.243$3E-01 

1/SEC 

0(40) 

0.2410SE-0) 

1/SEC 

0(41) 

O.277i4E-0l 

1/SEC 

0(42) 

0.233S1E-01 

1/$EC 

0(43) 

0.30800C-01 

1/SEC 

0(44) 

0.32143E-01 

l/SEC 

0(«$) 

0.33413E-01 

1/SEC 

0(44) 

e.34S4)C-0t 

l/SEC 

0(47) 

0.3$58(e-0t 

1/SEC 

0(48) 

0.34430E-01 

1/SEC 

0(43) 

0.37271E-01 

1/SEC 

0(50) 

0.37324E-01 

l/SEC 

G($n 

0.58«$3E'01 

1/SEC 

0(52) 

0.388$0E-D1 

l/SEC 

G(SS) 

0.33115E-0I 

1/SEC 

0(S4) 

0.53248E-01 

l/SEC 

ERAYS  •  F 

PRTALL 

•  F 

YO  •  1000.0  M 

DEPTH  • 

2000.0  H 

FREOC  «  250.0  HZ 

CHAX  • 

1500.8  M/SEC 

EREOC/CMAX  •  0.14447  CYCLES/M 

RATIO  • 

0.3333 

FRHAX  •  0.14465  CYCLES/M 

NFRS  a 

S 

DLTFR  ■  0.053330  CYCLES/M 

NOTEi  FRMIN  •  DLTFR 

RNGSTP 

■  1.0  M 

HRZRNO  ■  3000.0  M 

YR  a 

40.0  H 

NFR 

BETAO(DEG) 

DEPTH(M) 

TRVLT(SEC) 

BETAY(DEO) 

THETAT(RAO) 

THETAR(RAD) 

THETAY(RAD) 

MODULO  TWOPl 

TMETAY(RAO) 

1 

11.341 

134.18 

10.218510 

148.504 

0.1405120E*05 

0.4282557E*03 

0.1$42234E‘05 

4.0040532 

2 

23.140 

37.28 

5.108340 

154.443 

0.8025104E*04 

0.1254511E*04 

0.4748533E*04 

1.4024783 

3 

34.153 

1047.34 

3.40421$ 

34.157 

0.5350470E*04 

0. 1884747E«04 

0.3445702E>04 

.  3.4473828 

4 

$1,848 

478.77 

2.557817 

127.377 

0.4017810C«04 

8.2513023E*04 

0. 150478 7E*04 

3.1058334 

5 

73.434 

1522.23 

2.040408 

81.317 

0.3234735E*04 

0.3141278E*04 

0.3551484E*02 

1.243044$ 

total  CPU  TIME  •  0  MIN  ,  2(.(1  SEC 
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TABLE  3.18 

SIXTY-FIVE  DATA  PAIRS,  PIECEWISE  LINEAR  SOLUTION 


QAV  TRACIMO  USIMO  PIECEWISE  LINEAR  APPROXIMATION 
CACE  IS65 


ncata 

45 

NUHBEP 

DP  gradients 

•  44 

C(  1  ) 

-0.S*:56E-01 

I /SEC 

G(  2  1 

-0 . S4i40£-01 

l/SEC 

G(  S  ) 

-0  .  58471E-01 

1/SEC 

G(  C  1 

-0.!8*e8£-0l 

l/SEC 

C<  5) 

-0.38S1:E-01 

1/SEC 

G(  41 

-0 . 57«44e-01 

l/SEC 

G(  7) 

-0.S7:8«E-0I 

1/SEC 

G(  8  1 

-0.34435E-01 

l/SEC 

0(  1 ) 

-0  .  358«7E-01 

l/SEC 

G(  101 

-0.35073E-OI 

l/SEC 

C(  11  ) 

-0. SA 145E-0) 

l/SEC 

G(12) 

-0.1S174E-OI 

l/SEC 

GU  3  ) 

-0.3:iC!E-01 

l/SEC 

G(  14  1 

-0.!0»S5e-01 

1/SCC 

G(  15  ) 

-0.247!:E-01 

l/SEC 

G(  14  1 

-0.284386-01 

l/SEC 

G(  1  7  ) 

-0.27075E-01 

l/SEC 

G(  18  1 

-0.25447E-01 

l/SEC 

G(  !•  5 

-o.2<i5ee-oi 

l/SEC 

G(20  1 

-0  22410E-01 

l/SEC 

GC  1  ) 

-0.21007E-01 

l/SEC 

GCCl 

-0.145546-01 

l/SEC 

CCS  ) 

-0.  1  745<E-01 

l/SEC 

C(24  ) 

-0.  15412E-01 

l/SEC 

CCS) 

-0. i«i3:e-oi 

l/SEC 

G(24  1 

-0. I2317E-01 

l/SEC 

GC7) 

-O  .  10A73E-01 

l/SEC 

G(281 

-0.84032E-02 

)/SEC 

Gt:?) 

-0.471 JOE-02 

l/SEC 

G(301 

-0.480446-02 

I/ScC 

Gisi) 

-0.28884E-02 

l/SEC 

G(S2I 

-0.f4544£-0S 

l/SEC 

G(SS) 

0.9«34<e-03 

l/SEC 

0(34  1 

0. 2888 46-02 

1/SCC 

G(551 

0.<8044e-02 

l/SEC 

G(S4) 

0.47130E-02 

l/SEC 

6(371 

o.84C3:e-o: 

l/SEC 

C(!81 

0. 106736-01 

l/SEC 

G',9) 

o-irswE-oi 

l/SEC 

0(40) 

0.  14U2E-01 

l/SEC 

CUl  1 

0. 15112E-01 

l/SEC 

0(4:1 

0.17454E-01 

l/SEC 

Gt<5) 

0.  143S4E-01 

l/SEC 

0(44  1 

0.210076-01 

l/SEC 

G(«5  1 

0.2:410E-0I 

l/SEC 

C(44  > 

0.24I58E-0I 

l/SEC 

0(<.71 

0.254A7E-01 

wsec 

G(48  1 

0.270756-01 

l/SEC 

G(A11 

0.r8438E-0l 

l/SEC 

C(S0  1 

0.24;S36-OI 

l/SEC 

G(51  ) 

0.  S0«55e-01 

l/SEC 

0(52  1 

0.321036-01 

l/SEC 

0(5!  ) 

0.5J17AE-01 

l/SEC 

0(54) 

0.341456-01 

l/SEC 

Gt55) 

0.35073E-01 

l/SEC 

GI541 

0.558»7E-Cl 

l/SEC 

G(57) 

0,!44!5£-0l 

l/SEC 

G(S8) 

0.3?:84E-0l 

l/SEC 

0!5f  1 

0.378A4e-0l 

l/SEC 

0(40) 

0.383126-01 

l/SEC 

Gt«l  } 

0  .  38488E-0I 

1/SCC 

0(42) 

0.38«71E-0| 

l/SEC 

G(4S) 

ERavs 

0. 34140E-01 

F 

i/sec 

6(44}  • 

pbta 

0.5«:54E-0l 

LL  •  F 

l/SEC 

YO  »  1000 .0  M 

depth  •  2000.0  M 

PPEOC  •  250.0  M? 

PPEOC/CMAX  •  0. 14447 

CVCLES/W  rati 

0  •  0.4444 

FRMAX  •  0.14445  CVCLES/M 

DLYFR 

NOTE  r 

»  O.OSJSSO  CVCLES/M 

PRMJN  •  OLTFR 

RMGS 

TP  .  1 .0  M 

MRiRNG  •  3000.0  M 

CMAX  .  1500.0  M/SEC 

MFRS  •  5 

YR  .  A  0 . 0  M 


NPR 

B6TA0(DEC) 

DEPTHIH) 

TRVLTCSEC) 

BCTAVfOEC) 

thetatcrad) 

TMETAR(RAD) 

THETAY(RAD) 

MODULO  TWOPI 

TMETAV(RAD] 

1 

U  541 

134.17 

10.2185:8 

148  504 

0  .  1405123E*05 

0.4282557E*03 

0. 154:247E»05 

4.0340889 

2 

23.140 

37.27 

5  -  108444 

154  444 

0.80251 1BE*04 

0.1254511E»04 

0.47484076*04 

1 .4143483 

3 

34.  153 

1047.44 

3.4C4::i 

34  .  157 

0.5350474E‘04 

0.18847476*04 

0.34457126*04 

3.4744515 

4 

51 .848 

478 . 74 

2.5579:: 

127.478 

0.40i;8l8E*04 

0.:513025€*04 

0. 1504;45E.04 

S  .  1  1  3425  1 

5 

7«.4*4 

1522 . 3  ; 

2.04041: 

81.317 

0 . 3::48Q2E*04 

0.31412786*04 

0.«55:340E*0: 

1.2758237 

TCT 

AU  CPL  time 

»  0  MIN 

.  24.84  SEC 
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TABLE  3.19 

INCREASED  RANGE  STEP,  ODE  SOLVER  SOLUTION 


RAv  TRACING  USING  AKIHA  CUBIC  SPLINE  t  ODE  SOLVER 
CASE  :si7 

NDATA  •  17 


ERAVS  •  P 

PRTALL  •  P 

VO  •  1000.0  M 

DEPTH  •  2000.0  « 

FREOC  •  250.0  H2 

CMAX  >  1500.0  M/SEC 

PREOC/CMAX  ■  0.14*47  CVCLES/M 

RATIO  * 

FRMAX  •  0.U445  CVCLES/M 

NFRS  ■  5 

DLTFR  .  0.0S33SC  CVCLES/M 

NOTE:  PRHIN  •  DLTFR 

RNGGTP  •  S.O  M 

HR2RNG  •  SOOO.O  M 

VR  m  AO.O  M 

MFR 

6E7A0(DCG) 

DEPTm(N) 

TRVLT(SEC) 

BETAVfDEG) 

THETATIRAO) 

THETAR(RAD) 

THETAYIRAD) 

MODULO  TWOPI 

THETAYtRAD) 

1 

11.341 

134.04 

10.218584 

148.504 

0.  U05)S1E*05 

0.4282557E»03 

0. 1542304E*05 

4 .1214844 

2 

23.  140 

37.25 

5.  108  4  75 

154.444 

0.8025154E*04 

0.  1254StlE*04 

0.4748648E*04 

1.4570457 

3 

54  .  153 

1047 . 44 

3.404234 

34.157 

0.5350507E»04 

0. 1884747E*04 

0.3445740E*04 

3.7047770 

4 

51 .840 

478 . 74 

2.557834 

127. 478 

0.401783PE«04 

0.2513023E*04 

0 . 15048 1 7E«04 

3.1353500 

5 

74.446 

1522  .  34 

2 . 040625 

81.317 

0. 3r54822E*04 

0.314I278E*04 

0.4554S?'»E»C2 

1 .2460153 

TOTAL  CPU  TIME  •  10  HIM  .  5-53  SEC 


and  hundreds  or  thousandths  of  a  radian  for  modulo  27r  6^,  the 
reduction  in  CPU  time  is  following  a  nearly  linear  relation 
to  the  increase  in  range  step,  that  is,  CPU  time  was  reduced 
approximately  by  one-fifth. 

F.  EIGENRAYS 

Tables  3.20  and  3.21  show  the  eigenrays  found  by  the  two 
ray  acoustics  methods.  Figures  3.13  and  3.14  correspond  to 
these  two  tables.  The  figures  show  the  sound-speed  profile  to 
be  that  of  a  SOFAR  channel.  This  channel  profile  was  chosen 
because  it  produces  interesting  ray  patterns,  and  it  commonly 
occurs  in  nature.  The  minimum  sound-speed  occurs  at  1000 
meters  depth  and  a  local  minima  occurs  at  the  ocean  surface. 
The  tables  show  that  each  method  found  the  same  eigenrays  at 
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TABLE  3.20 

EIGENRAYS  FOUND  USING  THE 
PIECEWISE  LINEAR  APPROXIMATION  METHOD 


RAY  TRACINO  riEC£Wlt:E  LINEAR  APPROX  1 MA  t  I OM 

CAie  ICH-i 

NOATA  »  <  MUMPER  OP  GRAOIENtS  ’  5 

Gt  1)  «  O.UOOOE-OJ  1/CEC  Z)  *  -  0 . 1808  •»£- 0 1  I/SEC 

G(  3)  »  O.WOOOE-01  l/GEC 

ERAV$  9  r  PRTALL  »  P 

VO  »  100.0  M  DEPTH  •  lOOO.O  M  FREOC  »  256.0  MZ 

RWG5TP  ■  5.0  M  AMGGTP  »  l.O  DEC  HftZRMG  •  200.0  M 

VERROR  •  5 . 0  M 

eigemrays 

RAV  PETAOTDEG)  DEPTHfHJ  TRVLTfGEC)  BETAVttEG)  IMETaTIRAO)  THEIARIRAO) 

1  53.000  251.11  0.U7r*«  52.056  0 . 262 ;462e *05  0 . 1 6 724 • I E *  0 S 

2  54.000  245. 6^  0.165U2  53.P56  0 . 25  6  !56  ‘OS  0  •  I  6*42S2E  *  OS 

3  150.000  246.73  0.26  7040  24. ^SA  0 . 4  I  •/.644E*05  0  .  I  04  70  fSE  *  0  5 

TOTAL  CPU  TIME  ■  0  HIM  .  22.4?  GEC 


TABLE  3.21 

EIGENRAYS  FOUND  USING  THE 
ODE  SOLVER  METHOD 


RAV  TRACING  UIING  AVIMA  CUBIC  GPLINE  8  ODE  SOLVER 
CAGE  2CH4 

MCAfA  ■  4 

jOAv;  .  T  rpfAtt.  f  P 

VO  r  lOO.O  M  DEPTH  -  2000.0  M  PPEOr  •  250.0  HZ 

RNUIfP  *  5.0  M  AH7-,»r  »  1.0  pEG  MnZRMG  •  200.0  M 

VERPOO  f  5.0  M 

F  ir.FHRAvS 

OAV  PETACHDEG)  DEPIMfMJ  TOVLTIGEC)  PETAVfrEG)  IHETAf(RAD)  THETARIRAD) 

1  53.0cn  250  05  0.167065  52  «50  0 . 2 6242 1 AE *03  0 . I  6 7 24  *  I E • 0 3 

2  54  000  245.44  0  16-.«12  53.»52  0 . 25 '*04  55E  •  0  3  0  .  I  6  *42  !:E  *  0  3 

3  150.000  246.56  0.2*-*3'.J  2*. *71  0 . 4  j  ■»  1 56  3E  *03  0  - 1  04  70  •  SE  •  0  3 

total  CPU  TjMr  ,  hIM  ,  12. *4  'EC 


VR  «  250. 0  M 


THETAYIRAD)  HOrULO  TWQPI 
THETAVtRAO) 
0.*51725PE‘02  0.«24723« 

0  0*62Oj:e*o:  1.65572:4 

0.314447PE‘03  0.2077732 


VR  •  250.0  M 


THETAVCRAD1  MODULO  TWOPI 
THETAV(RAD) 
0.954«?ieE*02  I  .24A3187 

0  8**537  1E‘02  I  .  96*11  15 

O.SU7556E‘OS  0.5963779 
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launch  angles  of  53,  54  and  150  degrees-.  The  tables  highlight 
this  problem's  parameters  to  be: 

-  a  source  depth  yg  of  100  meters; 

-  a  receiver  depth  of  250  meters; 

-  and  a  horizontal  range  (HRZRNG)  to  the  receiver  of  200 
meters. 

As  expected,  the  ODE  solver  runs  into  range  limitations 
imposed  by  CPU  time  requirements.  The  200  meter  example 
presented  consumed  over  24  minutes  of  CPU  time  using  the  ODE 
solver  versus  the  21  seconds  used  by  the  piecewise  linear 
method . 

Tables  3.22  and  3.23  (along  with  their  corresponding 
Figures  3.15  and  3.16)  show  solutions  to  long-range  eigenray 
search  problems.  These  two  examples  continue  to  use  a  source 
depth  yg  of  100  meters  and  a  receiver  depth  y^  of  2  50  meters. 
Table  3.22  shows  a  horizontal  range  (HRZRNG)  of  10  kilometers 
with  a  depth  error  y^^ror  receiver  of  2  meters.  The 

angle  step  size  (angstp)  between  each  ray  is  0.1°.  Table  3.22 
also  shows  that  three  eigenrays  were  found  with  launch  angles 
of  51.2°,  90.9°  and  129.5°.  The  CPU  time  for  this  10km  run  is 
shown  to  be  8.5  minutes. 
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TABLE  3.22 

LONG  RANGE  EIGENRAYS,  PIECEWISE  LINEAR  SOLUTION 


RAV  TRACING  UGINO  PIECEWJSE  LINEAR  APPRCXIMAT ION 
cage  ICHC 

hdata  =  A  NUMBER  OP  C-RAOrENTS  *  J 

G«  n  =  O.UOOCE-Ol  l/CEC  G{  Z)  »  -0.ia88®£-0l  l/GEC 
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TABLE  3.23 

LONG  RANGE  EIGENRAYS,  PIECEWISE  LINEAR  SOLUTION 


RAV  TRACING  UGING  PieCENISE  LINEAR  APPROXIMATION 
CAGE  ICMA 

N04TA  s  <  NUMBER  OP  GRADIENTS  •  5 

G(  n  •  0.1«OCC£-01  l/GEC  C(  D  i  -0.l898«e-01  l/SEC 

C(  3  )  •  0. I  700CE  01  1/GEC 

ERAVS  *  T  PRTALL  *  P 

VO  «  100.0  M  DEPTH  *  rooo.o  M  FOEoc  «  rso.o  H? 

RNGGTP  >  B.O  M  ANOGTP  »  l.O  DEG  HR^RNG  •  50000.0  M 

VEPROR  »  15 . 0  M 

EIGEndavs 

RAV  BETAOirEG)  DEPTHIH)  TRVLnSECT  BETAvipEG)  TMPtaTIRAOT  TMETARIRAD) 

I  60.000  z:6.6^  52.6ri9;3  S».»17  0.8:36601E*05  0 . 5355 C *  IE ♦ 05 

Z  nr. 000  :51.00  36.660^**?  icr.535  0.560'»8'>«E*05  0.5iri05rE*05 

3  in. 00c  :-9.35  39.5r08ir  €C.8Z0  0.505?935E*05  0 . 65  7«06  OE  •  05 

6  i::.coo  :57.76  s**. 736555  irr.ier  0.5r6i507E*05  o.663««:5E‘05 

total  CPU  time  -  A  MIN  .  r7.51  SEC 


vR  s  rso. 


THETAVIRAD) 

0.6e5»l0?E*05 
O.r880U7E*O6 
0. 16717«5E‘05 
0  .  1801581E«05 


M 


MODULO  TWOPI 

THETAVIPADl 

6 . 1  3505:8 
2 . 7:8:732 
I .  oroso?® 


48 
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£IGEMRflY3:3  YERRCR:  2.3 


Table  3.23  shows  simulation  results  using  the  piecewise 
linear  method  with  a  horizontal  range  (HRZRNG)  of  50km. 
Parameters  that  affect  execution  time  such  as  angle  step  size, 
depths,  etc.,  match  the  simulation  presented  in  Table  3.21 
using  the  ODE  solver.  Table  3.21  shows  that  the  ODE  solver 
used  24  minutes  to  complete  a  200  meter  range  eigenvalue 
problem.  In  comparison.  Table  3.23  shows  that  the  piecewise 
linear  method  completed  the  50km  run  in  only  4  minutes. 

In  summary,  this  chapter  has  shown  that 

-  the  ODE  solver  performs  accurate  phase  calculations 
with  far  fewer  data  pair  samples  when  a  smooth 
sound-speed  profile  is  encountered; 

-  increasing  the  range  step  size  of  the  ODE  solver  to 
five  meters  reduces  the  CPU  time  cost  without 
signficantly  affecting  the  phase  solution  accuracies; 

-  the  piecewise  linear  method  can  run  long-range 
simulations  using  relatively  little  CPU  time;  and 

-  the  piecewise  linear  method  requires  many  sound-speed 
data  pair  samples  for  an  accurate  solution  in  a  medium 
with  a  smooth,  sinusoidal  profile. 
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IV.  CONCLUSIONS  AND  RECOMMENDATIONS 


The  transfer  function  of  the  inhomogeneous  ocean  based  on 
the  WKB  approximation  requires  solving  a  phase  integral.  Ray 
acoustics  theory  can  provide  solutions  to  this  phase  integral. 
The  simulations  performed  show  that  the  travel  times 
calculated  using  the  theories  of  ray  acoustics  can  be  used  to 
solve  the  phase  integral  avoiding  direct  numerical 
integration. 

Two  applications  of  ray  acoustics  produced  computer 
simulation  codes  which 

-  are  capable  of  solving  for  the  position,  travel  time 
and  phase  of  a  propagating  ray,  and 

-  have  very  different  advantages  and  costs. 

The  first  application  was  the  piecewise  linear  approximation. 
Sound-speed  versus  depth  data  pairs  sampled  from  the  ocean 
medium  were  connected  with  constant  gradient  linear  segments. 
Well-known,  closed  form  equations  form  the  mathematical  model 
for  sound  propagation.  The  solutions  are  low  cost  (in  terms 
of  CPU  time) ,  but  many  data  samples  are  required  for  accurate 
phase  solutions  for  arbitrary  sound-speed  profiles. 

In  contrast,  the  Akima  cubic  spline/ODE  solver  method  uses 
the  medium  samples  to  form  a  continuously  variable  sound-speed 
profile.  Accurate  phase  calculations  can  be  made  with  a 
minimum  of  sound-speed  versus  depth  data  pairs.  The 
disadvantage  in  using  the  ODE  solver  is  its  exorbitant  cost  in 
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terms  of  CPU  time.  Simulations  must  be  short  range  problems 
to  limit  computing  costs. 

The  ability  to  search  for  and  identify  eigenrays  was 
developed  for  each  of  the  two  ray  acoustics  theory 
applications.  Very  distinctive  characteristics  were  seen  for 
phase  calculations  in  terms  of  the  number  of  data  samples 
needed  to  assure  accurate  solutions.  The  identification  of 
eigenrays  is  an  easier  task  of  position  or  depth  computation. 
Here  the  distinctions  tend  to  disappear  with  the  two  methods 
arriving  at  the  same  solutions  for  eigenray  launch  angles. 
The  cost  or  CPU  time  required  continues  to  strongly  favor 
using  the  piecewise  linear  approach. 

These  findings  indicate  that  most  propagation  problems 
will  require  a  piecewise  linear  approach  for  computational 
efficiency.  Using  the  ODE  solver  would  severely  limit  the 
range  of  computer  simulations. 

Careful  sampling  of  the  ocean  media  is  required  to  obtain 
accurate  results  from  the  piecewise  linear  method.  Sufficient 
sound-speed  versus  depth  data  pairs  must  be  used  to  accurately 
represent  the  sound-speed  profile. 

This  thesis  developed  the  tools  to  quantify  the  strengths 
and  weaknesses  of  two  phase  computation  methods  in  a  variety 
of  media.  With  these  findings  in  mind,  future  work 
recommendations  are  to 

-  incorporate  each  phase  computation  technique  as  a 
module  in  the  larger  pulse  propagation  code;  and 
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-  run  pulse  propagation  simulations  to  compare  the 
received  pulse  shapes. 

If  this  comparison  shows  that  the  piecewise  linear 
approximation  produces  a  relatively  undistorted  received 
pulse,  it  is  an  efficient  solution  to  the  phase  computation 
problem. 
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